home *** CD-ROM | disk | FTP | other *** search
- #include <stdio.h>
- #include <math.h>
- #include "geom.h"
- #include "geomclass.h"
- #include "transform.h"
- #include "sphereP.h"
- #include "hpoint3.h"
-
- #ifndef FUDGE
- #define FUDGE 1e-6
- #endif
-
- void SphereMinMax(Sphere *sphere, HPoint3 *min, HPoint3 *max)
- {
- Geom *bbox;
- bbox = GeomBound((Geom *)sphere, TM_IDENTITY);
- BBoxMinMax((BBox *)bbox, min, max);
- GeomDelete(bbox);
- HPt3Normalize(min, min);
- HPt3Normalize(max, max);
- }
-
- void SphereCenter(Sphere *sphere, HPoint3 *center)
- {
- *center = sphere->center;
- }
-
- float SphereRadius(Sphere *sphere)
- {
- return sphere->radius;
- }
-
- void SphereSwitchSpace(Sphere *sphere, int space) {
- Transform T;
-
- sphere->space = space;
- TmScale(sphere->axis, sphere->radius, sphere->radius, sphere->radius);
- TmSpaceTranslate(T, sphere->center.x, sphere->center.y, sphere->center.z,
- sphere->space);
- GeomTransform((Geom*)sphere, T);
- }
-
- Sphere *SphereUnion3(Sphere *a, Sphere *b, Sphere *dest) {
- Sphere *sphere;
- int space;
- HPoint3 center, diff;
- float radius, mag_diff;
-
- if (a == NULL && b == NULL) return NULL;
- space = (a != NULL) ? a->space : b->space;
- if (dest == NULL)
- dest = (Sphere *)GeomCreate("sphere", CR_SPACE, space, CR_END);
-
- if (a == NULL || b == NULL) {
- if (a == NULL) {
- radius = b->radius;
- center = b->center;
- space = b->space;
- }
- else if (b == NULL) {
- radius = a->radius;
- center = a->center;
- space = a->space;
- }
- GeomSet((Geom*)dest, CR_RADIUS, radius, CR_CENTER, ¢er, CR_SPACE, space,
- CR_END);
- } else {
-
- if (a->space != b->space)
- OOGLError(1, "Uniting two spheres existing in different spaces.");
- space = a->space;
- if (space != TM_EUCLIDEAN)
- OOGLError(1, "SphereUnion3 currently only works reliably in\n%s",
- "Euclidean space.");
-
- HPt3Sub(&b->center, &a->center, &diff);
- Pt3Unit((Point3*)&diff);
- center.x = b->center.x + diff.x*b->radius;
- center.y = b->center.y + diff.y*b->radius;
- center.z = b->center.z + diff.z*b->radius;
- center.w = 1.0;
-
- GeomSet((Geom*)dest, CR_RADIUS, a->radius, CR_CENTER, &a->center, CR_END);
- SphereAddHPt3(dest, ¢er, TM_IDENTITY);
- }
-
- return dest;
- }
-
- void SphereEncompassBounds(Sphere *sphere, HPoint3 *points) {
- int i, j;
- float span, maxspan;
- HPoint3 *d1, *d2, center;
-
-
- d1 = d2 = &points[0];
- maxspan = 0.0;
- for (i = 0; i < 6; i++)
- for (j = i+1; j < 6; j++) {
- span = HPt3SpaceDistance(&points[i], &points[j], sphere->space);
- if (span > maxspan) {
- maxspan = span;
- d1 = &points[i];
- d2 = &points[j];
- }
- }
-
- /* Find the midpoint here - this will not work in non-Euclidean space */
- center.x = (d1->x + d2->x) / 2.0;
- center.y = (d1->y + d2->y) / 2.0;
- center.z = (d1->z + d2->z) / 2.0;
- center.w = 1.0;
-
- GeomSet((Geom *)sphere, CR_RADIUS, maxspan / 2.0, CR_CENTER, ¢er,
- CR_END);
-
- }
-
- int SphereAddHPt3(Sphere *sphere, HPoint3 *point, Transform T)
- {
- float radius, old_to_p, old_to_new;
- Point3 direction;
- HPoint3 center, newpoint;
-
- HPt3Transform(T, point, &newpoint);
- old_to_p = HPt3SpaceDistance(&newpoint, &sphere->center, sphere->space);
- if (old_to_p > sphere->radius) {
- radius = (sphere->radius + old_to_p) / 2.0;
- old_to_new = old_to_p - radius;
- center.x = sphere->center.x +
- (newpoint.x - sphere->center.x)*old_to_new / old_to_p;
- center.y = sphere->center.y +
- (newpoint.y - sphere->center.y)*old_to_new / old_to_p;
- center.z = sphere->center.z +
- (newpoint.z - sphere->center.z)*old_to_new / old_to_p;
- center.w = 1.0;
- GeomSet((Geom *)sphere, CR_RADIUS, radius, CR_CENTER, ¢er, CR_END);
- return 1;
- } else return 0;
- }
-
- int SphereAddHPt3N(Sphere *sphere, HPoint3 *point, int n, Transform T)
- {
- int i, ans = 0;
-
- for (i = 0; i < n; i++) ans |= SphereAddHPt3(sphere, &point[i], T);
- return ans;
- }
-
- void SphereEncompassHPt3N(Sphere *sphere, HPoint3 *point, int n,
- Transform T) {
- int i;
- HPoint3 spanPts[6];
-
- if (!n) return;
- for (i = 0; i < 6; i++) spanPts[i] = point[0];
- MaxDimensionalSpanN(spanPts, point, n);
- HPt3TransformN(T, spanPts, spanPts, 6);
- SphereEncompassBounds(sphere, spanPts);
- SphereAddHPt3N(sphere, point, n, T);
- }
-
- /* 0 has min x, 1 has max x, 2 has min y... */
- void MaxDimensionalSpan(HPoint3 *spanPts, HPoint3 *point)
- {
- if (point->x < spanPts[0].x) spanPts[0] = *point;
- else if (point->x > spanPts[1].x) spanPts[1] = *point;
- if (point->y < spanPts[2].y) spanPts[2] = *point;
- else if (point->y > spanPts[3].y) spanPts[3] = *point;
- if (point->z < spanPts[4].z) spanPts[4] = *point;
- else if (point->z > spanPts[5].z) spanPts[5] = *point;
- }
-
-
- void MaxDimensionalSpanN(HPoint3 *spanPts, HPoint3 *points, int n)
- {
- int i;
-
- for (i = 0; i < n; i++) {
- if (points[i].x < spanPts[0].x) spanPts[0] = points[i];
- else if (points[i].x > spanPts[1].x) spanPts[1] = points[i];
-
- if (points[i].y < spanPts[2].y) spanPts[2] = points[i];
- else if (points[i].y > spanPts[3].y) spanPts[3] = points[i];
-
- if (points[i].z < spanPts[4].z) spanPts[4] = points[i];
- else if (points[i].z > spanPts[5].z) spanPts[5] = points[i];
- }
- }
-